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Abstract 

The critical temperature and condensate fraction of a trapped interact- 
ing Bose gas are investigated when both atom-atom interaction and finite-size 
effects are taken into account. Canonical ensemble is used to obtain the equa- 
tions on the condensate fraction for the trapped interacting Bose gas near and 
below the critical temperature. In our approaches corrections due to atom- 
atom interaction and finite-size effects are obtained simultaneously for the 
critical temperature and condensate fraction of the system. Analytical high- 
order correction to the condensate fraction is given in this work. 
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I. INTRODUCTION 



The experimental realization of Bose-Einstein condensation (BEC) in the dilute alkali- 
metal atoms Q and more recently in the atomic hydrogen M has stimulated a new interest 
in the theoretical study of the inhomogeneous Bose gas. Thermodynamic properties such as 
critical temperature, condensate fraction have been investigated by several authors for the 
trapped Bose gas recently pi]. 

In the presence of an anisotropic harmonic potential of the 
form V ext (r) = m{uj 2 x x 2 + u 2 y 2 + u 2 z 2 )/2 the noninteracting model gives the critical tem- 
perature T c ° = ^ H? where Uho = (^x^y^z) 1 ^ 3 is the geometrical average of the 
oscillator frequencies, ((n) is Riemann ( function. In the large-N limit, the condensate frac- 
tion is £ = ^jr = 1 — • Finite-size effects [|J and interatomic interaction will give 
corrections to the thermodynamic properties of the system. The correction to the critical 
temperature due to finite-size effects has shown to be <5T C °/T C ° ~ — O^StoN" 1 / 3 /uho, where 
ZJ = {oj x + uj y + u z )/3 is the mean frequency. 

Although the atom clouds realized in the present experiments are very dilute, the effects 
due to interatomic interaction are important at low temperature. Researches show that 
atom-atom interaction will give leading corrections when N is larger than 10 5 for the alkalis. 
In fact the question of how two-body forces affect the thermodynamic properties of these 
systems have been the object of several theoretical investigations ||. Using a local density 
approximation Giorgini, Pitaevskii, and Stringari |||| obtained the shift of the critical tem- 
perature due to interatomic interaction: <5T C °/T C ° ~ — 1.33aA^ 1 / 6 /a/ lo , where ato = Jfijnuoho 
is the harmonic oscillator length. Recently Monte Carlo simulation is also used to inves- 
tigate the thermodynamic properties of the trapped interacting Bose gases. 

In this work, we investigate the critical temperature and the condensate fraction of 
the system when both atom-atom interaction and finite-size effects are taken into account. 
Canonical ensemble is used to obtain the analytical high-order correction due to interatomic 
interaction. Especially, finite-size effects are obtained simultaneously. To obtain the shift of 
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the critical temperature, we give the analytic description of the condensation fraction near 
the critical temperature. 

The paper is planed as follows. In Sec. II we outline the canonical ensemble. In Sec. Ill 
we give the condensate fraction of the system near the critical temperature. The shift of the 
critical temperature agrees with the usual results In Sec. IV we obtain the analytical 

high-order correction to the condensate fraction due to atom-atom interaction. The lowest- 
order correction agrees with the well-established results [||,|]||. Sec.V the correction to 
the condensate fraction beyond mean field theory is given. 



II. PARTITION FUNCTION OF THE TRAPPED INTERACTING BOSE GASES 

Canonical ensemble and saddle-point method have been used to investigate the ther- 
modynamic properties of the interacting uniform Bose gases ||. In this work canonical 
ensemble is used to discuss the trapped interacting Bose gases. The partition function of N 
trapped interacting bosons is given by 

Q(N) = S n exp {~f3E n ) = ^{n x ,n v ,n z } exp [-(3 \^n x ,n y ^n x n y n^n x n y n z + E int }] > W 

where N nxnynz and e nxnynz are the occupation numbers and energy level of the state 
{n x , n y , n z } respectively. E int is the interaction energy of the system. 
From (D 

Q{N)= ££ o=0 {exp [-(3 (E + Ei„t)] Qo (N - N )} , (2) 

where Q (N - iV ) = Z' {nx ,n y ,n z } ex P \~P (^n Xl n y ,n z N nxnyriz e nxnynz )] stands for the partition 
function of a fictitious system of iV — Nq trapped noninteracting atoms. Eq is the energy of 
the condensate. For convenience, we have separated out the ground state n x = n y = n z = 
from the states {n x ,n y ,n z } ^ 0, so that we first carry out the primed summation over all 
{n x ,n y ,n z } ({n x ,n y ,n z } ^ 0) with a fixed value (N — N ) of the partial sum ^\ nxnynz } an d 
then carry out a summation over all possible values of N , ie. N = to N = N. Assume 
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A (N - N ) is the free energy of the fictitious system. A (N - N ) = -± lnQ (iV - N ). 
From (|) 

Q(N)= S^ o=0 exp [-/3 (E + E mt ) - (3A (N - N )\ . (3) 

The sum S]y _q in (|3]) maybe replaced by the largest term in the sum, for the error omitted 
in doing so will be statistically negligible. With this approximation we do not investigate 
the fluctuations of the system. However, this approximation is reasonable because we can 
give the occupation number of the energy level s nxTlynz , which agrees with the widely used 
mean occupation number in the frame of grand-canonical ensemble. Neglecting all terms 
but the largest in (^), the number of bosons in the condensate can be obtained as: 

~ P dN {E ° + Emt) ~ P dN MN ~ No) = °- (4) 
The calculations of the free energy Aq(N — No) of the fictitious noninteracting Bose gas 
is nontrivial because there is a requirement that the number of the particles is N — No in the 
summation of the partition function. Using the saddle-point method of integration developed 
by Darwin and Fowler |T(| it is straightforward to find that — (3-t^Ao(N — N ) = lnz in 
which zq is the fugacity of the fictitious Bose gas. In addition, the fugacity z is determined 
by 

N - N o= K M — r- R — 1 i — r - ( 5 ) 

eXD ft£n-n,,n, Zr, 1 



We can easily understand (|^) in terms of Bose-Einstein distribution of the trapped non- 



interacting Bose gas. Using the relation — /3-^A (N — iVo) = lnzo> (1) becomes 

- (3^- (E + E int ) + lnz = 0. (6) 

01\q 

Equations (||) and (|) will be used to discuss the corrections due to atom-atom interaction 
and finite-size effects. Once we know the interaction energy of the system, it is easy to obtain 
the correction due to atom-atom interaction. In addition, finite-size effects are separated off 
in (|5]). Omitting interactions between atoms, we can obtain \nz = /3eooo from (|5]). From ([5]) 
the number of the thermal atoms Nt is given by 



N T = N- N = E' ^ -T—. (7) 



exp 



l3{Sn x n y n z ~ ^000, 
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This is the exact conclusion in the frame of grand-canonical ensemble for the trapped 
noninteracting Bose gas. This shows the equivalence between the grand-canonical ensemble 
and canonical ensemble in the calculations of the condensate fraction. In addition, it shows 
that neglecting all terms but the largest in (§) is reasonable. In the absence of interaction, 
from (0), the condensate fraction is given by || 

£ = 1 - f— ) 3 - 3 ^ C(2) (— ) 2 N-^ 3 (8) 

The third term in is the finite-size correction to the ideal Bose gas in the thermody- 
namic limit. In short, (||) and ([]) account for both interatomic interaction and finite-size 
effects. We will use (Q) and (^|) to calculate the critical temperature and condensate fraction 
of the system. 

III. CRITICAL TEMPERATURE AND CONDENSATE FRACTION NEAR T c 

A. interaction energy of the system near the critical temperature 

The parameter expressing the importance of the interatomic interaction compared to the 
kinetic energy is JN 1 oc ^H. ; where a is the scattering length between bosons. Near the 
critical temperature << 1. This means that E^ in >> E int . In this case we can use the 



method of pseudopotentials developed by Huang and Yang ]TT[ to calculate the interaction 
energy of the system. 

In the method of pseudopotentials, the actual Hamiltonian of the system is replaced by 
an effective Hamiltonian, such that the ground state and the low-lying energy levels of the 
system are given equally well by the new Hamiltonian. With hard-sphere approximation, the 
boundary conditions between atoms are replaced by the pseudopotential operator g8^ {f)4-r 
in the effective Hamiltonian, where g = 4?ra ^ 2 . The effective Hamiltonian of the system may 
be taken to be 



h 2 

H = -— (V? + • • • + V%) + Z? =1 V ext (r t ) + S Ki ^, (9) 

where Uij = gH^jS^ 3 ' (r*j — r}) For a trapped dilute Bose gas near the critical tem- 

perature the pseudopotentials can be regarded as perturbation terms. The energy levels to 
the first order in a may be obtained through the usual perturbation theory. 

In the general case of a system containing an arbitrary number N of bosons, the unper- 
turbed normalized wave function of the system is given by 

(n m i--.\ 1/2 

* N= [ N\' ) (10) 
where the sum is taken over all permutations of the different suffixes a±, c*2, • • •, cxn and 
the numbers N a . show how many of these suffixes have the same value aij (with T,N ai = N). 
From the first-order perturbation theory the interaction energy E int of the system is given 
by 

E int =< ^ N \^i<j^ij\^ N > • (11) 

Thus, the operators (d/drij) will be operating on a set of functions which are well 
behaved for all values of r^-; accordingly, these operators may be replaced by the unit 
operators. Basing on the investigation of the various permutations |J of the single particle 
states the interaction energy takes the form, 

E mt = 9 ^~ No) f <^(f)d 3 f + 2gN f <^(r> T (f) d 3 f + g f n 2 T (r)d 3 r. (12) 
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To obtain fll2"D we used the fact that ny(f) = S^ o ^0f (f) is the density distribution of 



the thermal atoms. Near the critical temperature, the density distribution of the condensate 
is n if) = N (pl (f) = iVo { m ^ L y' /2 ex P\~T i^^x 2 + u y y 2 + uj z z 2 )}. The density distribution 
of the thermal atoms can be obtained using Bose-Einstein distribution and the semiclas- 
sical approximation of the energy level of a single particle. Near the critical temperature 
^■j— oc iV 1//3 . In the available traps N ranges from a few thousand to several millions, thus 
fiuJho « ksT®, ie. the energy level of the thermal atoms can be approximated as con- 
tinuous. In addition, this means that the semiclassical approximation for the normal gas 
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is expected to work well on a wide range of temperatures 0. In terms of Bose-Einstein 
distribution n T (r) = J d 3 p(2irh)~ 3 [exp (f3e (p,r)) — l] -1 , where e(p,r) = J^- + V^t (r) is 
the semiclassical energy in phase space and (3 = 1/ksT. Density distribution of the ther- 
mal atoms in the coordinate space would be ^(r) = Ay 3 g 3 /2 (exp[— V ext (r)/(kBT)]), where 
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> V 2 

2irh I (m/csT) is the thermal wavelength. (z) belongs to the class of func- 
tions g a (z) = T^ =l z n /n a . By integrating the semiclassical approximation nj-(r) over space 
one obtain the condensate fraction £ = 1 — (^5-) , which does not account for finite-size 
effects. However, in this paper the semiclassical approximation riT{r) is used to calculate 
the interaction energy of the system, ie. used to calculate the high-order correction due 
to interatomic interaction. Thus the adoption of the semiclassical approximation will omit 
only higher-order modification, which is of the order of the multiplication of the corrections 
due to the interatomic interaction and finite-size effects. Because the finite-size effects have 
been separated off in (||), we can give finite-size effects although semiclassical approximation 
is used to calculate the correction due to interatomic interaction. 

Using the geometrical average of the oscillator frequencies Uho, ut{t) = 



At 3 <?3/2 (exp [— muj\ r 2 1'lksTX)- Introducing the width of the normal gas Rt = y 2&bT /ma> 2 
and defining the rescaled variable through r = Rt^, we can obtain tit (r) = 
(exp [— r 2 ]). Near the critical temperature ksT >> hcuho, the classical Boltzmann 
distribution ti c i (r) oc Ay exp [— r 2 ] gives a zero-order approximation for the density distri- 
bution of the thermal cloud. Thus we take the trial function of the form 

n T (?) = Vi^t 3 ex P ~V2^ 2 (13) 

where 771 and 772 are two dimensionless variational parameters. Another reason for the 
adoption of the Gaussian trial function lies in the fact that there is a factor exp [— f 2 ] in 
tit (r). rji and r] 2 are determined by the requirement that J \n^f (r) — (r) \d 3 r obtain 
the minimum. Results of numerical calculations show that 7/1 = 2.587, 772 = 1.146 and 
J \n^f (f) —riT (r) \d 3 f/ J nx (r) d 3 f= 0.0017. So the Gaussian distribution n!f (r) agrees very 
well with nT(r)- We should note that the Gaussian distribution n£ (r) is different from the 



classical Boltzmann distribution because of the variational parameters rji and 772- 
Combining (0) and (|X3|) we have 

E mt = 1_ - 6 K + (2 - ^L) r)ig\^ 3 N + M^iV, (14) 
2 I v 2,71(1^0 ) 



where we have used the fact that a^ >> and /c^T » htu near the critical temperature. 

B. shift of the critical temperature and condensate fraction near T c 

We introduce a scaling parameter which accounts for the role of the two-body repulsive 
interaction. This parameter 9 is fixed by the ratio between griT (r = 0,T C °) and the critical 



temperature k B T® . 9 = grtT ( r J^f ) _ jV 1 / 6 . The scaling parameter 6* can also be 

kb^c C(3) 

written in the form = 0.65r7 5/2 , where 77 = /i (T = 0) /k B T® = C{3)1/2 2 152/5 (^ 1/6 ^;) 2/5 is 
also an important scaling parameter accounting for the role of interatomic interaction P,|12" . 
Obviously # stands for high-order correction due to interatomic interaction, compared to the 
parameter 77. From (|^), (^), and (]I4|) we can obtain the equation on the condensate fraction. 

'-^wf^+S*^ 1 ' 3 ' (15) 

_c(^ f _ (2 __y wl/2+ln ^ 0i (16) 

where we have used the reduced temperature t = T/T®. From (|1|) and (|16|) the condensate 
fraction is given by 

1-t 3 - 4.51-2-^/6^/2 _ ^2^-1/3 

5 ~ 1 + 1.16-2-iV 2 / 3 * 2 ' 1 ' 

a ho 

The third term in the numerator of ( |TTD represents the correction due to the interaction 
between atoms, while the last term in the numerator accounts for the correction due to 
the finite-size effects. By setting £ = in ( |l~7|) one can obtain the shift of the critical 
temperature, 

S = S J™ + 5 Il^l = -1.50— JVV« - ^ 1/3 - (18) 

T c ° TO TO a ho 2C(3) 2/3 ^o 1 ' 
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The first term in fll8|) is the shift due to interatomic interaction. It agrees with the 
results based on the local density approximation Because of the denominator in (p~7|) . 
close to T c ° the condensate fraction will increase slowly with decreasing the temperature. 
Especially, it shows that there is only high-order correction of the scaling parameter rj due 
to interatomic interaction. We should note that ( JP7| ) holds only when << 1. The second 
contribution in ( |l8l) gives exactly the usual results due to the finite-size effects. 



IV. CONDENSATE FRACTION OF THE SYSTEM BELOW T c 

A. lowest-order correction due to interatomic interaction and finite-size effects 

Below the critical temperature >> 1, it. E int » E^ in . In this case we can use 
the well-known Thomas-Fermi approximation. With the Thomas-Fermi approximation the 
density profile of the condensate is n (f) = — ^ _ y ext (f)). The normalization 
condition on no(r) provides the relation between chemical potential and number of particles 
in the condensate: fi = f 1 ^° a ) ^ ■ From Gross-Pitaevskii (GP) equation JOJ and the 
well-known Virial theorem the energy of the condensate turns out to be Eq = (5/7) fiNo 
||. Omitting atom-atom interactions in the normal gas and the interaction between the 
condensate and normal gas, from @ and (|) 

N = N - N T = N - K x ,n y , nz TTT (19) 

exp \p [e nxnynz - /ij 

From (ID it is easy to obtain the following result 



£ = i _ f _ C(2) _£ 3a<(2) 2 1/3 , s 



where the last term in (|20| ) is the correction due to finite-size effects [[§]. Using the relation 
n/hsT = ?7^ 2//5 /t, the third term in (^) gives exactly the usual lowest-order modification 
due to atom-atom interactions ||,|6|,|8[] . 

C(3) ft 2c C(3) 2/3 
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This proves the equivalence between canonical ensemble and grand-canonical ensemble, 
even in the presence of interatomic interactions. We will give high-order modification in the 
following. 



B. high-order modification due to interatomic interactions 



Equation ( |21| ) is obtained by omitting the interaction between thermal atoms and inter- 
action between thermal atoms and the condensate. When these two sorts of interactions are 
considered 
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En + Ei, 



int 



7' 



fiN + 2 g / no (r) nx if) d r + g n T (r ) d f 



(22) 



The calculations of (|22|) are straightforward. The results is given by 



E + E, 



int 



7 



+ 
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35 X z T k B T 2 3 /2 



(23) 



where = J 2«/ mw^ is the radius of the condensate, i? 3 / A| = oir\ where 
ct RA = 0.294 (a/a ho ) 3/5 . From (||) 



0«iT (^o + 3 



From (^), (|]), and (^) we can obtain the condensate fraction below the critical temper 
ature, 



1 - 0.89r) o/z t 



5/2.1/2 



1.01r] o/ H 



5/2.1/2 



(24) 



e = i - 1 3 - 



cm 

C(3) 



(l - O.S^H 1 / 2 ) ^ + l.Olry 5 /^ 1 / 2 



3a7C(2) 
2^ C (3) 2/3 



t 2 N -i/s (25) 



In (|25 ) the terms comprises ifl 2 represents high-order correction to the condensate frac- 
tion due to interatomic interaction. Omitting the terms comprises 7/ 5//2 , we can obtain the 
lowest-order correction fl2~l"|). Essentially, (^1|) and fl25|) are transcendental equations on the 
condensate fraction. We can easily give the numerical solution of the condensate fraction. 



In Fig.l and Fig. 2 we used the experimental parameter by Ensher et al. | I4]| , where the 
cloud consists of 4 x 10 4 atoms at the transition and aja^ = 5.4 x 10~ 3 . According to 
Fig. 2 the high-order correction (p5|) agrees well with the condensate fraction (|17]) in the 
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region << 1. From (|T7|) and (|25|) , close to T c ° the condensate fraction increase quite 



slowly with decreasing the temperature. However, the lowest-order correction (£]]) is rather 
different from the dashed line near the critical temperature. The difference between ( |P7D 
and fl21| ) lies in the fact that there is only high-order correction due to interatomic interac- 
tion near the critical temperature, while there is only lowest-order correction in Fig. 2 
demonstrates the validity of (^) over the whole range of temperature below T c . 

In the preceding calculations we neglects the role of the interaction between the conden- 
sate and the thermal atoms on the density distributions of the condensate and the thermal 
atoms. In the frame of Hartree-Fork model ||J|, the densities of the condensate and the 
thermal component are given as follows. 

n (f) = ^ye Xt {r)-2gn T {f) e ^ _ ^ (f) _ ^ (26) 

U T (f) = A r 3 (7 3/2 [ e -(^ t (^+29(noM+n T (r))- M )/fc s Tj _ ^ 

Using (H), (H), and (pBf), (|27|) and the usual iterative procedure we can obtain the 
numerical conclusion of the condensate fraction. The numerical conclusion is illustrated in 
Fig.3. 

V. CONDENSATE FRACTION BEYOND MEAN-FIELD THEORY 

In the preceding calculations of the condensate fraction, GP equation is used to obtain 
the energy E of the condensate. GP equation is expected to be valid if the system is dilute, 
ie. n\a\ 3 << 1. Obviously, correction to GP equation would modify the condensate fraction 
of the system. The first correction to the mean-field approximation have been investigated 



by Timmermans, Tommasini, and Huang [15] and by Braaten and Nieto [16]. For large iVo, 



using the local density approximation, the density distribution of the condensate is given by 

4 m 3/2 

3ir 2 h 2 



V-V ext (r) 4m 3 / 2 /2 
no (r) = — -2 \p ~ Vext (r)J , (28) 



with jj, given by 
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hu ho /15N a\ 



(^)>y«). 

The parameter a 3 n (0) can be directly expressed in terms of the relevant parameters of 
the system. a 3 n (0) = {^^S^ = 0.0078r? 6 £ 2 / 5 . Using the relation E = 5N fx/7 
and (0), (M) we can obtain the correction to the condensate fraction below T c 



C(2),3 
C(3)* 



1 - 0.89r? 5/2 t 1/2 + 0.16^ 



t 



3u7C(2) 
2uU (3) 2/3 



The term O.Wi] 3 ^ 1 ^ 5 is the correction due to the modification of the GP equation. It gives 
higher-order correction of the scaling parameter rj, compared to the high-order correction 
due to interatomic interactions. 

Another correction to the condensate fraction beyond mean-field approximation is given 
by the quantum depletion of the condensate. There are atoms which do not occupy the 
condensate at zero temperature because of quantum depletion effects. We can use local 
density approximation to write the density of atoms out of the condensate. Timmermans, 
Tommasini, and Huang gives n out (f) = (8/3) [n (f) a 3 /^] 1 ^ 2 . Integration of n out (r) yields 
the result 



iout = J n out if) d 3 r 



a 3 n(0). 



(31) 



Combining ( j30|) and flST]) 



e=l-t 3 



C(2),3 
C(3)* 



o.sgrf/Y/ 2 + o.i6^ 1/5 ) ^ + lmrf'H 

3cJC (2) 



1/2 



0.11r/ 3 C 1/5 - 



(32) 



2u ho C (3) 2/3 

The term — 0.1 lr] 3 ^ 1 / 5 is the correction due to quantum depletion. 

Compared to the finite-size effects, the corrections to the mean-field approximation have 
a different dependence on the parameter N and a/a^. The corrections beyond mean field 
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theory become larger than finite-size effects when N is larger than about 10 6 . In Fig. 4 
we give numerical conclusion of (|32|). The corrections beyond mean field theory become 
important when the temperature is much lower than the critical temperature. 

VI. CONCLUSION 

To conclude, the canonical ensemble is used to obtain the equations on the condensate 
fraction (|5]) and which account for interaction effects. Finite-size effects can be obtained 
simultaneously because these effects are comprised in (||). @ and @ are used to investigate 
the condensate fraction near and below the critical temperature. From the condensate 
fraction near the critical temperature, we obtain the shift of the critical temperature due 
to atom-atom interaction and finite-size effects simultaneously. In addition, the analytical 
high-order correction due to interatomic interaction is obtained in this work. 
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FIGURES 

FIG. 1. Theoretical prediction for the condensate fraction vs T/T® for N = 4 x 10 4 rubidium 
atoms in a trap with a/a/j G = 5.4 x 10 -3 . The dashed line shows the noninteracting model in the 
thermodynamical limit. The dotted line is obtained from the numerical solution of which 
accounts for the lowest-order correction due to interatomic interaction. The solid line is obtained 
from (p5[) which accounts for high-order correction. 



FIG. 2. The upper solid line shows the lowest-order correction due to interatomic inter- 
action. The lower solid line represents the high-order correction (pq). The dashed line shows the 
condensate fraction near T c , which is obtained from (|i~7l). Equation (|i~7|) holds when << 1. In 



the region — — << 1, the lower solid line agrees well with the dashed line. This demonstrates the 



Npa 

a ho 

validity of (p^) over the whole range of temperature T <T C , even close to T®. 



FIG. 3. The solid line is obtained from pq). The dashed line is a self-consistent numerical 
solution of (||), @, (|2"2|), and (pq 1 ), fl27|) using iterative method. 

FIG. 4. The solid line is obtained from (p5f) which accounts for high-order correction due to 



interatomic interaction. The dotted line is obtained from (32), in which the corrections beyond 
mean field theory are comprised. The corrections beyond mean-field theory become important 
when the temperature is much lower than the critical temperature (shown in the inset). 
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